#pragma ident "$Id: $"

/**
 * @file HarrisPriesterDrag.cpp
 * The HarrisPriester class computes the Harris-Priester atmosphere model.
 */

//============================================================================
//
//  This file is part of GPSTk, the GPS Toolkit.
//
//  The GPSTk is free software; you can redistribute it and/or modify
//  it under the terms of the GNU Lesser General Public License as published
//  by the Free Software Foundation; either version 2.1 of the License, or
//  any later version.
//
//  The GPSTk is distributed in the hope that it will be useful,
//  but WITHOUT ANY WARRANTY; without even the implied warranty of
//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//  GNU Lesser General Public License for more details.
//
//  You should have received a copy of the GNU Lesser General Public
//  License along with GPSTk; if not, write to the Free Software Foundation,
//  Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
//
//  Wei Yan - Chinese Academy of Sciences . 2009, 2010
//
//============================================================================


#include "HarrisPriesterDrag.hpp"
#include "IERS.hpp"
#include "ReferenceFrames.hpp"
#include "StringUtils.hpp"
#include <map>

namespace gpstk
{
   using namespace std;


      // Upper height limit [km]
   const double HarrisPriesterDrag::upper_limit = 2000.0;         
      
      // Lower height limit [km]
   const double HarrisPriesterDrag::lower_limit = 100.0;  

      // Right ascension lag [rad]
   const double HarrisPriesterDrag::ra_lag      = 0.523599;          

      // Harris-Priester atmospheric density coefficients
      // it's modified from JAT by Richard C.
   const double HarrisPriesterDrag::hpcoef[10][177] = {
   
   // ALTITUDE              MINIMUM DENSITY       MAXIMUM DENSITY

   {   // F107 = +65.000000
      +1.1000000e+002,      +7.8000000e+001,      +7.8000000e+001,
      +1.2000000e+002,      +2.4900000e+001,      +2.4900000e+001,
      +1.3000000e+002,      +8.9780000e+000,      +9.3310000e+000,
      +1.4000000e+002,      +4.0690000e+000,      +4.2120000e+000,
      +1.5000000e+002,      +2.0860000e+000,      +2.1680000e+000,
      +1.6000000e+002,      +1.1460000e+000,      +1.2360000e+000,
      +1.7000000e+002,      +6.6160000e-001,      +7.5580000e-001,
      +1.8000000e+002,      +4.0160000e-001,      +4.8850000e-001,
      +1.9000000e+002,      +2.5300000e-001,      +3.2740000e-001,
      +2.0000000e+002,      +1.6280000e-001,      +2.2840000e-001,
      +2.1000000e+002,      +1.0760000e-001,      +1.6340000e-001,
      +2.2000000e+002,      +7.2870000e-002,      +1.1920000e-001,
      +2.3000000e+002,      +5.0380000e-002,      +8.8510000e-002,
      +2.4000000e+002,      +3.5490000e-002,      +6.6660000e-002,
      +2.5000000e+002,      +2.5410000e-002,      +5.0830000e-002,
      +2.6000000e+002,      +1.8460000e-002,      +3.9190000e-002,
      +2.7000000e+002,      +1.3580000e-002,      +3.0500000e-002,
      +2.8000000e+002,      +1.0100000e-002,      +2.3940000e-002,
      +2.9000000e+002,      +7.5880000e-003,      +1.8940000e-002,
      +3.0000000e+002,      +5.7190000e-003,      +1.5100000e-002,
      +3.2000000e+002,      +3.3050000e-003,      +9.8860000e-003,
      +3.4000000e+002,      +1.9530000e-003,      +6.6080000e-003,
      +3.6000000e+002,      +1.1750000e-003,      +4.4940000e-003,
      +3.8000000e+002,      +7.1670000e-004,      +3.1000000e-003,
      +4.0000000e+002,      +4.4280000e-004,      +2.1630000e-003,
      +4.2000000e+002,      +2.7790000e-004,      +1.5260000e-003,
      +4.4000000e+002,      +1.7600000e-004,      +1.0850000e-003,
      +4.6000000e+002,      +1.1280000e-004,      +7.7670000e-004,
      +4.8000000e+002,      +7.3460000e-005,      +5.5990000e-004,
      +5.0000000e+002,      +4.8660000e-005,      +4.0610000e-004,
      +5.2000000e+002,      +3.2910000e-005,      +2.9630000e-004,
      +5.4000000e+002,      +2.2790000e-005,      +2.1740000e-004,
      +5.6000000e+002,      +1.6220000e-005,      +1.6050000e-004,
      +5.8000000e+002,      +1.1880000e-005,      +1.1920000e-004,
      +6.0000000e+002,      +8.9780000e-006,      +8.9100000e-005,
      +6.2000000e+002,      +6.9870000e-006,      +6.7080000e-005,
      +6.4000000e+002,      +5.5930000e-006,      +5.0900000e-005,
      +6.6000000e+002,      +4.5890000e-006,      +3.8960000e-005,
      +6.8000000e+002,      +3.8460000e-006,      +3.0110000e-005,
      +7.0000000e+002,      +3.2810000e-006,      +2.3510000e-005,
      +7.2000000e+002,      +2.8380000e-006,      +1.8570000e-005,
      +7.4000000e+002,      +2.4820000e-006,      +1.4840000e-005,
      +7.6000000e+002,      +2.1900000e-006,      +1.2020000e-005,
      +7.8000000e+002,      +1.9440000e-006,      +9.8570000e-006,
      +8.0000000e+002,      +1.7360000e-006,      +8.1930000e-006,
      +8.5000000e+002,      +1.1800000e-006,      +6.2000000e-006,
      +9.0000000e+002,      +8.7000000e-007,      +4.4000000e-006,
      +9.5000000e+002,      +6.6000000e-007,      +3.3000000e-006,
      +1.0000000e+003,      +4.8000000e-007,      +2.7000000e-006,
      +1.1000000e+003,      +3.0000000e-007,      +1.7500000e-006,
      +1.2000000e+003,      +1.8500000e-007,      +1.2000000e-006,
      +1.3000000e+003,      +1.1300000e-007,      +8.5000000e-007,
      +1.4000000e+003,      +7.3000000e-008,      +6.2000000e-007,
      +1.5000000e+003,      +5.2000000e-008,      +4.7500000e-007,
      +1.6000000e+003,      +3.7000000e-008,      +3.6500000e-007,
      +1.7000000e+003,      +2.5500000e-008,      +3.0000000e-007,
      +1.8000000e+003,      +1.8200000e-008,      +2.2000000e-007,
      +1.9000000e+003,      +1.3000000e-008,      +1.8000000e-007,
      +2.0000000e+003,      +1.0000000e-008,      +1.4500000e-007,
   },
   {   // F107 = +75.0000000
      +1.1000000e+002,      +7.8000000e+001,      +7.8000000e+001,
      +1.2000000e+002,      +2.4900000e+001,      +2.4900000e+001,
      +1.3000000e+002,      +8.9130000e+000,      +9.2630000e+000,
      +1.4000000e+002,      +4.0510000e+000,      +4.1970000e+000,
      +1.5000000e+002,      +2.0920000e+000,      +2.1760000e+000,
      +1.6000000e+002,      +1.1600000e+000,      +1.2490000e+000,
      +1.7000000e+002,      +6.7810000e-001,      +7.6970000e-001,
      +1.8000000e+002,      +4.1570000e-001,      +5.0160000e-001,
      +1.9000000e+002,      +2.6530000e-001,      +3.3910000e-001,
      +2.0000000e+002,      +1.7250000e-001,      +2.3810000e-001,
      +2.1000000e+002,      +1.1510000e-001,      +1.7160000e-001,
      +2.2000000e+002,      +7.8660000e-002,      +1.2620000e-001,
      +2.3000000e+002,      +5.4850000e-002,      +9.4350000e-002,
      +2.4000000e+002,      +3.8940000e-002,      +7.1550000e-002,
      +2.5000000e+002,      +2.8100000e-002,      +5.4920000e-002,
      +2.6000000e+002,      +2.0560000e-002,      +4.2610000e-002,
      +2.7000000e+002,      +1.5230000e-002,      +3.3370000e-002,
      +2.8000000e+002,      +1.1410000e-002,      +2.6350000e-002,
      +2.9000000e+002,      +8.6310000e-003,      +2.0960000e-002,
      +3.0000000e+002,      +6.5540000e-003,      +1.6820000e-002,
      +3.2000000e+002,      +3.8400000e-003,      +1.1130000e-002,
      +3.4000000e+002,      +2.3010000e-003,      +7.5140000e-003,
      +3.6000000e+002,      +1.4030000e-003,      +5.1600000e-003,
      +3.8000000e+002,      +8.6840000e-004,      +3.5950000e-003,
      +4.0000000e+002,      +5.4400000e-004,      +2.5340000e-003,
      +4.2000000e+002,      +3.4190000e-004,      +1.7980000e-003,
      +4.4000000e+002,      +2.1890000e-004,      +1.2910000e-003,
      +4.6000000e+002,      +1.4180000e-004,      +9.3340000e-004,
      +4.8000000e+002,      +9.3040000e-005,      +6.7950000e-004,
      +5.0000000e+002,      +6.1940000e-005,      +4.9760000e-004,
      +5.2000000e+002,      +4.1950000e-005,      +3.6640000e-004,
      +5.4000000e+002,      +2.8990000e-005,      +2.7130000e-004,
      +5.6000000e+002,      +2.0500000e-005,      +2.0190000e-004,
      +5.8000000e+002,      +1.4870000e-005,      +1.5110000e-004,
      +6.0000000e+002,      +1.1090000e-005,      +1.1370000e-004,
      +6.2000000e+002,      +8.5000000e-006,      +8.6140000e-005,
      +6.4000000e+002,      +6.6990000e-006,      +6.5670000e-005,
      +6.6000000e+002,      +5.4160000e-006,      +5.0430000e-005,
      +6.8000000e+002,      +4.4790000e-006,      +3.9030000e-005,
      +7.0000000e+002,      +3.7770000e-006,      +3.0470000e-005,
      +7.2000000e+002,      +3.2380000e-006,      +2.4020000e-005,
      +7.4000000e+002,      +2.8120000e-006,      +1.9130000e-005,
      +7.6000000e+002,      +2.4680000e-006,      +1.5400000e-005,
      +7.8000000e+002,      +2.1840000e-006,      +1.2540000e-005,
      +8.0000000e+002,      +1.9450000e-006,      +1.0330000e-005,
      +8.5000000e+002,      +1.4000000e-006,      +7.3000000e-006,
      +9.0000000e+002,      +1.0000000e-006,      +5.2000000e-006,
      +9.5000000e+002,      +7.8000000e-007,      +3.8000000e-006,
      +1.0000000e+003,      +5.8000000e-007,      +3.1000000e-006,
      +1.1000000e+003,      +3.5500000e-007,      +2.0500000e-006,
      +1.2000000e+003,      +2.3000000e-007,      +1.4000000e-006,
      +1.3000000e+003,      +1.4500000e-007,      +1.0000000e-006,
      +1.4000000e+003,      +9.1000000e-008,      +7.2000000e-007,
      +1.5000000e+003,      +6.7000000e-008,      +5.4500000e-007,
      +1.6000000e+003,      +4.5000000e-008,      +4.1500000e-007,
      +1.7000000e+003,      +3.3000000e-008,      +3.4000000e-007,
      +1.8000000e+003,      +2.2500000e-008,      +2.6500000e-007,
      +1.9000000e+003,      +1.6000000e-008,      +2.1000000e-007,
      +2.0000000e+003,      +1.2800000e-008,      +1.7000000e-007,
   },
   {   // F107 = +100.0000000
      +1.1000000e+002,      +7.8000000e+001,      +7.8000000e+001,
      +1.2000000e+002,      +2.4900000e+001,      +2.4900000e+001,
      +1.3000000e+002,      +8.7400000e+000,      +9.0910000e+000,
      +1.4000000e+002,      +4.0020000e+000,      +4.1570000e+000,
      +1.5000000e+002,      +2.1030000e+000,      +2.1910000e+000,
      +1.6000000e+002,      +1.1950000e+000,      +1.2830000e+000,
      +1.7000000e+002,      +7.1860000e-001,      +8.0690000e-001,
      +1.8000000e+002,      +4.5150000e-001,      +5.3530000e-001,
      +1.9000000e+002,      +2.9570000e-001,      +3.6950000e-001,
      +2.0000000e+002,      +1.9800000e-001,      +2.6380000e-001,
      +2.1000000e+002,      +1.3540000e-001,      +1.9380000e-001,
      +2.2000000e+002,      +9.4590000e-002,      +1.4520000e-001,
      +2.3000000e+002,      +6.7370000e-002,      +1.1050000e-001,
      +2.4000000e+002,      +4.8800000e-002,      +8.5270000e-002,
      +2.5000000e+002,      +3.5890000e-002,      +6.6570000e-002,
      +2.6000000e+002,      +2.6750000e-002,      +5.2500000e-002,
      +2.7000000e+002,      +2.0170000e-002,      +4.1770000e-002,
      +2.8000000e+002,      +1.5380000e-002,      +3.3490000e-002,
      +2.9000000e+002,      +1.1830000e-002,      +2.7050000e-002,
      +3.0000000e+002,      +9.1170000e-003,      +2.2020000e-002,
      +3.2000000e+002,      +5.5220000e-003,      +1.4960000e-002,
      +3.4000000e+002,      +3.4210000e-003,      +1.0360000e-002,
      +3.6000000e+002,      +2.1580000e-003,      +7.2980000e-003,
      +3.8000000e+002,      +1.3820000e-003,      +5.2100000e-003,
      +4.0000000e+002,      +8.9520000e-004,      +3.7630000e-003,
      +4.2000000e+002,      +5.8630000e-004,      +2.7440000e-003,
      +4.4000000e+002,      +3.8750000e-004,      +2.0190000e-003,
      +4.6000000e+002,      +2.5840000e-004,      +1.4970000e-003,
      +4.8000000e+002,      +1.7390000e-004,      +1.1170000e-003,
      +5.0000000e+002,      +1.1820000e-004,      +8.3820000e-004,
      +5.2000000e+002,      +8.1230000e-005,      +6.3240000e-004,
      +5.4000000e+002,      +5.6500000e-005,      +4.7940000e-004,
      +5.6000000e+002,      +3.9870000e-005,      +3.6510000e-004,
      +5.8000000e+002,      +2.8590000e-005,      +2.7930000e-004,
      +6.0000000e+002,      +2.0880000e-005,      +2.1460000e-004,
      +6.2000000e+002,      +1.5560000e-005,      +1.6560000e-004,
      +6.4000000e+002,      +1.1850000e-005,      +1.2840000e-004,
      +6.6000000e+002,      +9.2300000e-006,      +1.0000000e-004,
      +6.8000000e+002,      +7.3510000e-006,      +7.8300000e-005,
      +7.0000000e+002,      +5.9820000e-006,      +6.1620000e-005,
      +7.2000000e+002,      +4.9650000e-006,      +4.8790000e-005,
      +7.4000000e+002,      +4.1950000e-006,      +3.8860000e-005,
      +7.6000000e+002,      +3.6000000e-006,      +3.1170000e-005,
      +7.8000000e+002,      +3.1300000e-006,      +2.5180000e-005,
      +8.0000000e+002,      +2.7510000e-006,      +2.0500000e-005,
      +8.5000000e+002,      +2.0000000e-006,      +1.3000000e-005,
      +9.0000000e+002,      +1.5000000e-006,      +9.0000000e-006,
      +9.5000000e+002,      +1.2000000e-006,      +6.2000000e-006,
      +1.0000000e+003,      +9.4000000e-007,      +5.0000000e-006,
      +1.1000000e+003,      +5.8000000e-007,      +3.1000000e-006,
      +1.2000000e+003,      +3.7500000e-007,      +2.1000000e-006,
      +1.3000000e+003,      +2.5000000e-007,      +1.5000000e-006,
      +1.4000000e+003,      +1.7000000e-007,      +1.1000000e-006,
      +1.5000000e+003,      +1.1800000e-007,      +7.7000000e-007,
      +1.6000000e+003,      +8.2000000e-008,      +6.0000000e-007,
      +1.7000000e+003,      +6.0000000e-008,      +5.0000000e-007,
      +1.8000000e+003,      +4.3000000e-008,      +3.9000000e-007,
      +1.9000000e+003,      +3.2000000e-008,      +3.1000000e-007,
      +2.0000000e+003,      +2.4000000e-008,      +2.5000000e-007,
   },
   {   // F107 = +125.0000000
      +1.1000000e+002,      +7.8000000e+001,      +7.8000000e+001,
      +1.2000000e+002,      +2.4900000e+001,      +2.4900000e+001,
      +1.3000000e+002,      +8.5670000e+000,      +8.9200000e+000,
      +1.4000000e+002,      +3.9520000e+000,      +4.1130000e+000,
      +1.5000000e+002,      +2.1110000e+000,      +2.2020000e+000,
      +1.6000000e+002,      +1.2260000e+000,      +1.3140000e+000,
      +1.7000000e+002,      +7.5720000e-001,      +8.4110000e-001,
      +1.8000000e+002,      +4.8640000e-001,      +5.6720000e-001,
      +1.9000000e+002,      +3.2550000e-001,      +3.9890000e-001,
      +2.0000000e+002,      +2.2430000e-001,      +2.8890000e-001,
      +2.1000000e+002,      +1.5670000e-001,      +2.1590000e-001,
      +2.2000000e+002,      +1.1180000e-001,      +1.6440000e-001,
      +2.3000000e+002,      +8.1150000e-002,      +1.2720000e-001,
      +2.4000000e+002,      +5.9870000e-002,      +9.9780000e-002,
      +2.5000000e+002,      +4.4790000e-002,      +7.9130000e-002,
      +2.6000000e+002,      +3.3940000e-002,      +6.3360000e-002,
      +2.7000000e+002,      +2.6000000e-002,      +5.1160000e-002,
      +2.8000000e+002,      +2.0110000e-002,      +4.1630000e-002,
      +2.9000000e+002,      +1.5700000e-002,      +3.4090000e-002,
      +3.0000000e+002,      +1.2310000e-002,      +2.8090000e-002,
      +3.2000000e+002,      +7.6750000e-003,      +1.9530000e-002,
      +3.4000000e+002,      +4.8940000e-003,      +1.3850000e-002,
      +3.6000000e+002,      +3.1780000e-003,      +9.9760000e-003,
      +3.8000000e+002,      +2.0940000e-003,      +7.2820000e-003,
      +4.0000000e+002,      +1.3970000e-003,      +5.3760000e-003,
      +4.2000000e+002,      +9.6250000e-004,      +4.0330000e-003,
      +4.4000000e+002,      +6.5520000e-004,      +3.0340000e-003,
      +4.6000000e+002,      +4.4970000e-004,      +2.2990000e-003,
      +4.8000000e+002,      +3.1090000e-004,      +1.7540000e-003,
      +5.0000000e+002,      +2.1660000e-004,      +1.3460000e-003,
      +5.2000000e+002,      +1.5200000e-004,      +1.0380000e-003,
      +5.4000000e+002,      +1.0760000e-004,      +8.0430000e-004,
      +5.6000000e+002,      +7.6850000e-005,      +6.2600000e-004,
      +5.8000000e+002,      +5.5450000e-005,      +4.8920000e-004,
      +6.0000000e+002,      +4.0470000e-005,      +3.8380000e-004,
      +6.2000000e+002,      +2.9920000e-005,      +3.0220000e-004,
      +6.4000000e+002,      +2.2440000e-005,      +2.3880000e-004,
      +6.6000000e+002,      +1.7100000e-005,      +1.8940000e-004,
      +6.8000000e+002,      +1.3250000e-005,      +1.5070000e-004,
      +7.0000000e+002,      +1.0450000e-005,      +1.2040000e-004,
      +7.2000000e+002,      +8.3970000e-006,      +9.6580000e-005,
      +7.4000000e+002,      +6.8670000e-006,      +7.7770000e-005,
      +7.6000000e+002,      +5.7130000e-006,      +6.2890000e-005,
      +7.8000000e+002,      +4.8300000e-006,      +5.1090000e-005,
      +8.0000000e+002,      +4.1430000e-006,      +4.1700000e-005,
      +8.5000000e+002,      +3.1500000e-006,      +2.7000000e-005,
      +9.0000000e+002,      +2.2500000e-006,      +1.8000000e-005,
      +9.5000000e+002,      +1.7200000e-006,      +1.2000000e-005,
      +1.0000000e+003,      +1.4000000e-006,      +8.8000000e-006,
      +1.1000000e+003,      +8.5000000e-007,      +5.3000000e-006,
      +1.2000000e+003,      +5.4000000e-007,      +3.2000000e-006,
      +1.3000000e+003,      +3.8000000e-007,      +2.2500000e-006,
      +1.4000000e+003,      +2.8000000e-007,      +1.5000000e-006,
      +1.5000000e+003,      +1.9500000e-007,      +1.0300000e-006,
      +1.6000000e+003,      +1.3500000e-007,      +8.0000000e-007,
      +1.7000000e+003,      +1.0200000e-007,      +6.5000000e-007,
      +1.8000000e+003,      +7.4000000e-008,      +5.2000000e-007,
      +1.9000000e+003,      +5.6000000e-008,      +4.2500000e-007,
      +2.0000000e+003,      +4.4000000e-008,      +3.4000000e-007,
   },
   {   // F107 = +150.0000000
      +1.1000000e+002,      +7.8000000e+001,      +7.8000000e+001,
      +1.2000000e+002,      +2.4900000e+001,      +2.4900000e+001,
      +1.3000000e+002,      +8.3770000e+000,      +8.7100000e+000,
      +1.4000000e+002,      +3.8990000e+000,      +4.0590000e+000,
      +1.5000000e+002,      +2.1220000e+000,      +2.2150000e+000,
      +1.6000000e+002,      +1.2630000e+000,      +1.3440000e+000,
      +1.7000000e+002,      +8.0080000e-001,      +8.7580000e-001,
      +1.8000000e+002,      +5.2830000e-001,      +6.0100000e-001,
      +1.9000000e+002,      +3.6170000e-001,      +4.2970000e-001,
      +2.0000000e+002,      +2.5570000e-001,      +3.1620000e-001,
      +2.1000000e+002,      +1.8390000e-001,      +2.3960000e-001,
      +2.2000000e+002,      +1.3410000e-001,      +1.8530000e-001,
      +2.3000000e+002,      +9.9490000e-002,      +1.4550000e-001,
      +2.4000000e+002,      +7.4880000e-002,      +1.1570000e-001,
      +2.5000000e+002,      +5.7090000e-002,      +9.3080000e-002,
      +2.6000000e+002,      +4.4030000e-002,      +7.5550000e-002,
      +2.7000000e+002,      +3.4300000e-002,      +6.1820000e-002,
      +2.8000000e+002,      +2.6970000e-002,      +5.0950000e-002,
      +2.9000000e+002,      +2.1390000e-002,      +4.2260000e-002,
      +3.0000000e+002,      +1.7080000e-002,      +3.5260000e-002,
      +3.2000000e+002,      +1.0990000e-002,      +2.5110000e-002,
      +3.4000000e+002,      +7.2140000e-003,      +1.8190000e-002,
      +3.6000000e+002,      +4.8240000e-003,      +1.3370000e-002,
      +3.8000000e+002,      +3.2740000e-003,      +9.9550000e-003,
      +4.0000000e+002,      +2.2490000e-003,      +7.4920000e-003,
      +4.2000000e+002,      +1.5580000e-003,      +5.6840000e-003,
      +4.4000000e+002,      +1.0910000e-003,      +4.3550000e-003,
      +4.6000000e+002,      +7.7010000e-004,      +3.3620000e-003,
      +4.8000000e+002,      +5.4740000e-004,      +2.6120000e-003,
      +5.0000000e+002,      +3.9160000e-004,      +2.0420000e-003,
      +5.2000000e+002,      +2.8190000e-004,      +1.6050000e-003,
      +5.4000000e+002,      +2.0420000e-004,      +1.2670000e-003,
      +5.6000000e+002,      +1.4880000e-004,      +1.0050000e-003,
      +5.8000000e+002,      +1.0920000e-004,      +7.9970000e-004,
      +6.0000000e+002,      +8.0700000e-005,      +6.3900000e-004,
      +6.2000000e+002,      +6.0120000e-005,      +5.1230000e-004,
      +6.4000000e+002,      +4.5190000e-005,      +4.1210000e-004,
      +6.6000000e+002,      +3.4300000e-005,      +3.3250000e-004,
      +6.8000000e+002,      +2.6320000e-005,      +2.6910000e-004,
      +7.0000000e+002,      +2.0430000e-005,      +2.1850000e-004,
      +7.2000000e+002,      +1.6070000e-005,      +1.7790000e-004,
      +7.4000000e+002,      +1.2810000e-005,      +1.4520000e-004,
      +7.6000000e+002,      +1.0360000e-005,      +1.1900000e-004,
      +7.8000000e+002,      +8.4960000e-006,      +9.7760000e-005,
      +8.0000000e+002,      +7.0690000e-006,      +8.0590000e-005,
      +8.5000000e+002,      +4.8000000e-006,      +5.5000000e-005,
      +9.0000000e+002,      +3.3000000e-006,      +3.7000000e-005,
      +9.5000000e+002,      +2.4500000e-006,      +2.4000000e-005,
      +1.0000000e+003,      +1.9000000e-006,      +1.7000000e-005,
      +1.1000000e+003,      +1.1800000e-006,      +8.7000000e-006,
      +1.2000000e+003,      +7.5000000e-007,      +4.8000000e-006,
      +1.3000000e+003,      +5.3000000e-007,      +3.2000000e-006,
      +1.4000000e+003,      +4.1000000e-007,      +2.0000000e-006,
      +1.5000000e+003,      +2.9000000e-007,      +1.3500000e-006,
      +1.6000000e+003,      +2.0000000e-007,      +9.5000000e-007,
      +1.7000000e+003,      +1.6000000e-007,      +7.7000000e-007,
      +1.8000000e+003,      +1.2000000e-007,      +6.3000000e-007,
      +1.9000000e+003,      +9.6000000e-008,      +5.2000000e-007,
      +2.0000000e+003,      +7.3000000e-008,      +4.4000000e-007,
   },
   {   // F107 = +175.0000000
      +1.1000000e+002,      +7.8000000e+001,      +7.8000000e+001,
      +1.2000000e+002,      +2.4900000e+001,      +2.4900000e+001,
      +1.3000000e+002,      +8.2390000e+000,      +8.5770000e+000,
      +1.4000000e+002,      +3.8580000e+000,      +4.0210000e+000,
      +1.5000000e+002,      +2.1230000e+000,      +2.2180000e+000,
      +1.6000000e+002,      +1.2830000e+000,      +1.3630000e+000,
      +1.7000000e+002,      +8.2460000e-001,      +8.9970000e-001,
      +1.8000000e+002,      +5.5380000e-001,      +6.2510000e-001,
      +1.9000000e+002,      +3.8480000e-001,      +4.5140000e-001,
      +2.0000000e+002,      +2.7600000e-001,      +3.3620000e-001,
      +2.1000000e+002,      +2.0190000e-001,      +2.5690000e-001,
      +2.2000000e+002,      +1.4950000e-001,      +2.0080000e-001,
      +2.3000000e+002,      +1.1240000e-001,      +1.5930000e-001,
      +2.4000000e+002,      +8.5790000e-002,      +1.2800000e-001,
      +2.5000000e+002,      +6.6260000e-002,      +1.0400000e-001,
      +2.6000000e+002,      +5.1740000e-002,      +8.5220000e-002,
      +2.7000000e+002,      +4.0790000e-002,      +7.0410000e-002,
      +2.8000000e+002,      +3.2440000e-002,      +5.8580000e-002,
      +2.9000000e+002,      +2.6000000e-002,      +4.9040000e-002,
      +3.0000000e+002,      +2.0870000e-002,      +4.1310000e-002,
      +3.2000000e+002,      +1.3680000e-002,      +2.9910000e-002,
      +3.4000000e+002,      +9.1570000e-003,      +2.2010000e-002,
      +3.6000000e+002,      +6.2380000e-003,      +1.6430000e-002,
      +3.8000000e+002,      +4.3130000e-003,      +1.2410000e-002,
      +4.0000000e+002,      +3.0190000e-003,      +9.4750000e-003,
      +4.2000000e+002,      +2.1340000e-003,      +7.2930000e-003,
      +4.4000000e+002,      +1.5230000e-003,      +5.6650000e-003,
      +4.6000000e+002,      +1.0960000e-003,      +4.4330000e-003,
      +4.8000000e+002,      +7.9340000e-004,      +3.4910000e-003,
      +5.0000000e+002,      +5.7820000e-004,      +2.7660000e-003,
      +5.2000000e+002,      +4.2370000e-004,      +2.2030000e-003,
      +5.4000000e+002,      +3.1220000e-004,      +1.7630000e-003,
      +5.6000000e+002,      +2.3130000e-004,      +1.4170000e-003,
      +5.8000000e+002,      +1.7230000e-004,      +1.1430000e-003,
      +6.0000000e+002,      +1.2900000e-004,      +9.2580000e-004,
      +6.2000000e+002,      +9.7240000e-005,      +7.5230000e-004,
      +6.4000000e+002,      +7.3740000e-005,      +6.1330000e-004,
      +6.6000000e+002,      +5.6310000e-005,      +5.0140000e-004,
      +6.8000000e+002,      +4.3340000e-005,      +4.1110000e-004,
      +7.0000000e+002,      +3.3630000e-005,      +3.3800000e-004,
      +7.2000000e+002,      +2.6340000e-005,      +2.7870000e-004,
      +7.4000000e+002,      +2.0840000e-005,      +2.3030000e-004,
      +7.6000000e+002,      +1.6660000e-005,      +1.9080000e-004,
      +7.8000000e+002,      +1.3470000e-005,      +1.5850000e-004,
      +8.0000000e+002,      +1.1030000e-005,      +1.3200000e-004,
      +8.5000000e+002,      +7.8000000e-006,      +9.2000000e-005,
      +9.0000000e+002,      +5.1000000e-006,      +6.3000000e-005,
      +9.5000000e+002,      +3.4500000e-006,      +4.1000000e-005,
      +1.0000000e+003,      +2.6000000e-006,      +2.8500000e-005,
      +1.1000000e+003,      +1.6000000e-006,      +1.4000000e-005,
      +1.2000000e+003,      +1.0000000e-006,      +7.4000000e-006,
      +1.3000000e+003,      +7.0000000e-007,      +4.3000000e-006,
      +1.4000000e+003,      +5.4300000e-007,      +2.6000000e-006,
      +1.5000000e+003,      +4.0000000e-007,      +1.7000000e-006,
      +1.6000000e+003,      +2.9000000e-007,      +1.2000000e-006,
      +1.7000000e+003,      +2.3000000e-007,      +9.4000000e-007,
      +1.8000000e+003,      +1.7500000e-007,      +7.4000000e-007,
      +1.9000000e+003,      +1.4000000e-007,      +6.0000000e-007,
      +2.0000000e+003,      +1.0900000e-007,      +5.2000000e-007,
   },
   {   // F107 = +200.0000000
      +1.1000000e+002,      +7.8000000e+001,      +7.8000000e+001,
      +1.2000000e+002,      +2.4900000e+001,      +2.4900000e+001,
      +1.3000000e+002,      +8.0680000e+000,      +8.3810000e+000,
      +1.4000000e+002,      +3.8110000e+000,      +3.9640000e+000,
      +1.5000000e+002,      +2.1300000e+000,      +2.2220000e+000,
      +1.6000000e+002,      +1.3130000e+000,      +1.3860000e+000,
      +1.7000000e+002,      +8.5680000e-001,      +9.2780000e-001,
      +1.8000000e+002,      +5.9370000e-001,      +6.5310000e-001,
      +1.9000000e+002,      +4.2380000e-001,      +4.7720000e-001,
      +2.0000000e+002,      +3.1000000e-001,      +3.5990000e-001,
      +2.1000000e+002,      +2.2970000e-001,      +2.7780000e-001,
      +2.2000000e+002,      +1.7320000e-001,      +2.1940000e-001,
      +2.3000000e+002,      +1.3260000e-001,      +1.7590000e-001,
      +2.4000000e+002,      +1.0290000e-001,      +1.4280000e-001,
      +2.5000000e+002,      +8.0730000e-002,      +1.1710000e-001,
      +2.6000000e+002,      +6.3980000e-002,      +9.6890000e-002,
      +2.7000000e+002,      +5.1170000e-002,      +8.0800000e-002,
      +2.8000000e+002,      +4.1240000e-002,      +6.7840000e-002,
      +2.9000000e+002,      +3.3490000e-002,      +5.7350000e-002,
      +3.0000000e+002,      +2.7350000e-002,      +4.8850000e-002,
      +3.2000000e+002,      +1.8380000e-002,      +3.5950000e-002,
      +3.4000000e+002,      +1.2600000e-002,      +2.6870000e-002,
      +3.6000000e+002,      +8.7780000e-003,      +2.0350000e-002,
      +3.8000000e+002,      +6.2050000e-003,      +1.5600000e-002,
      +4.0000000e+002,      +4.4410000e-003,      +1.2070000e-002,
      +4.2000000e+002,      +3.2100000e-003,      +9.4170000e-003,
      +4.4000000e+002,      +2.3430000e-003,      +7.4110000e-003,
      +4.6000000e+002,      +1.7230000e-003,      +5.8740000e-003,
      +4.8000000e+002,      +1.2760000e-003,      +4.6860000e-003,
      +5.0000000e+002,      +9.5060000e-004,      +3.7600000e-003,
      +5.2000000e+002,      +7.1220000e-004,      +3.0330000e-003,
      +5.4000000e+002,      +5.3620000e-004,      +2.4580000e-003,
      +5.6000000e+002,      +4.0570000e-004,      +2.0010000e-003,
      +5.8000000e+002,      +3.0830000e-004,      +1.6350000e-003,
      +6.0000000e+002,      +2.3540000e-004,      +1.3410000e-003,
      +6.2000000e+002,      +1.8050000e-004,      +1.1040000e-003,
      +6.4000000e+002,      +1.3900000e-004,      +9.1140000e-004,
      +6.6000000e+002,      +1.0760000e-004,      +7.5470000e-004,
      +6.8000000e+002,      +8.3720000e-005,      +6.2670000e-004,
      +7.0000000e+002,      +6.5490000e-005,      +5.2180000e-004,
      +7.2000000e+002,      +5.1530000e-005,      +4.3550000e-004,
      +7.4000000e+002,      +4.0800000e-005,      +3.6430000e-004,
      +7.6000000e+002,      +3.2530000e-005,      +3.0550000e-004,
      +7.8000000e+002,      +2.6130000e-005,      +2.5670000e-004,
      +8.0000000e+002,      +2.1160000e-005,      +2.1620000e-004,
      +8.5000000e+002,      +1.3000000e-005,      +1.4000000e-004,
      +9.0000000e+002,      +7.9000000e-006,      +9.4000000e-005,
      +9.5000000e+002,      +5.0000000e-006,      +6.4000000e-005,
      +1.0000000e+003,      +3.5000000e-006,      +4.4000000e-005,
      +1.1000000e+003,      +2.1000000e-006,      +2.2000000e-005,
      +1.2000000e+003,      +1.3000000e-006,      +1.1000000e-005,
      +1.3000000e+003,      +9.0000000e-007,      +6.2000000e-006,
      +1.4000000e+003,      +6.8000000e-007,      +3.6000000e-006,
      +1.5000000e+003,      +5.0000000e-007,      +2.2500000e-006,
      +1.6000000e+003,      +3.7000000e-007,      +1.5500000e-006,
      +1.7000000e+003,      +3.0000000e-007,      +1.2000000e-006,
      +1.8000000e+003,      +2.3200000e-007,      +9.0000000e-007,
      +1.9000000e+003,      +1.8200000e-007,      +7.3000000e-007,
      +2.0000000e+003,      +1.4500000e-007,      +6.0000000e-007,
   },
   {   // F107 = +225.0000000
      +1.1000000e+002,      +7.8000000e+001,      +7.8000000e+001,
      +1.2000000e+002,      +2.4900000e+001,      +2.4900000e+001,
      +1.3000000e+002,      +7.9560000e+000,      +8.2620000e+000,
      +1.4000000e+002,      +3.7770000e+000,      +3.9290000e+000,
      +1.5000000e+002,      +2.1280000e+000,      +2.2210000e+000,
      +1.6000000e+002,      +1.3240000e+000,      +1.3970000e+000,
      +1.7000000e+002,      +8.7820000e-001,      +9.4480000e-001,
      +1.8000000e+002,      +6.1190000e-001,      +6.7050000e-001,
      +1.9000000e+002,      +4.3800000e-001,      +4.9380000e-001,
      +2.0000000e+002,      +3.2300000e-001,      +3.7510000e-001,
      +2.1000000e+002,      +2.4410000e-001,      +2.9130000e-001,
      +2.2000000e+002,      +1.8680000e-001,      +2.3170000e-001,
      +2.3000000e+002,      +1.4450000e-001,      +1.8700000e-001,
      +2.4000000e+002,      +1.1330000e-001,      +1.5280000e-001,
      +2.5000000e+002,      +8.9760000e-002,      +1.2620000e-001,
      +2.6000000e+002,      +7.1820000e-002,      +1.0510000e-001,
      +2.7000000e+002,      +5.7950000e-002,      +8.8200000e-002,
      +2.8000000e+002,      +4.7130000e-002,      +7.4520000e-002,
      +2.9000000e+002,      +3.8500000e-002,      +6.3470000e-002,
      +3.0000000e+002,      +3.1540000e-002,      +5.4390000e-002,
      +3.2000000e+002,      +2.1520000e-002,      +4.0480000e-002,
      +3.4000000e+002,      +1.4960000e-002,      +3.0590000e-002,
      +3.6000000e+002,      +1.0570000e-002,      +2.3420000e-002,
      +3.8000000e+002,      +7.5690000e-003,      +1.8130000e-002,
      +4.0000000e+002,      +5.4870000e-003,      +1.4170000e-002,
      +4.2000000e+002,      +4.0140000e-003,      +1.1160000e-002,
      +4.4000000e+002,      +2.9670000e-003,      +8.8630000e-003,
      +4.6000000e+002,      +2.2100000e-003,      +7.0900000e-003,
      +4.8000000e+002,      +1.6570000e-003,      +5.7070000e-003,
      +5.0000000e+002,      +1.2510000e-003,      +4.6200000e-003,
      +5.2000000e+002,      +9.4890000e-004,      +3.7590000e-003,
      +5.4000000e+002,      +7.2360000e-004,      +3.0740000e-003,
      +5.6000000e+002,      +5.5430000e-004,      +2.5240000e-003,
      +5.8000000e+002,      +4.2640000e-004,      +2.0810000e-003,
      +6.0000000e+002,      +3.2940000e-004,      +1.7220000e-003,
      +6.2000000e+002,      +2.5550000e-004,      +1.4300000e-003,
      +6.4000000e+002,      +1.9890000e-004,      +1.1910000e-003,
      +6.6000000e+002,      +1.5550000e-004,      +9.9500000e-004,
      +6.8000000e+002,      +1.2210000e-004,      +8.3350000e-004,
      +7.0000000e+002,      +9.6310000e-005,      +7.0000000e-004,
      +7.2000000e+002,      +7.6300000e-005,      +5.8930000e-004,
      +7.4000000e+002,      +6.0740000e-005,      +4.9730000e-004,
      +7.6000000e+002,      +4.8620000e-005,      +4.2060000e-004,
      +7.8000000e+002,      +3.9130000e-005,      +3.5640000e-004,
      +8.0000000e+002,      +3.1690000e-005,      +3.0270000e-004,
      +8.5000000e+002,      +2.0000000e-005,      +2.0000000e-004,
      +9.0000000e+002,      +1.2000000e-005,      +1.3500000e-004,
      +9.5000000e+002,      +7.7000000e-006,      +9.3000000e-005,
      +1.0000000e+003,      +5.1000000e-006,      +6.8000000e-005,
      +1.1000000e+003,      +2.8000000e-006,      +3.4000000e-005,
      +1.2000000e+003,      +1.6800000e-006,      +1.7000000e-005,
      +1.3000000e+003,      +1.1000000e-006,      +9.4000000e-006,
      +1.4000000e+003,      +8.1000000e-007,      +5.2000000e-006,
      +1.5000000e+003,      +6.0000000e-007,      +3.2000000e-006,
      +1.6000000e+003,      +4.5000000e-007,      +2.2000000e-006,
      +1.7000000e+003,      +3.6000000e-007,      +1.6000000e-006,
      +1.8000000e+003,      +2.8000000e-007,      +1.2000000e-006,
      +1.9000000e+003,      +2.2500000e-007,      +9.1000000e-007,
      +2.0000000e+003,      +1.7800000e-007,      +7.2000000e-007,
   },
   {   // F107 = +250.0000000
      +1.1000000e+002,      +7.8000000e+001,      +7.8000000e+001,
      +1.2000000e+002,      +2.4900000e+001,      +2.4900000e+001,
      +1.3000000e+002,      +7.8320000e+000,      +8.1210000e+000,
      +1.4000000e+002,      +3.7420000e+000,      +3.8850000e+000,
      +1.5000000e+002,      +2.1300000e+000,      +2.2200000e+000,
      +1.6000000e+002,      +1.3400000e+000,      +1.4100000e+000,
      +1.7000000e+002,      +9.0040000e-001,      +9.6240000e-001,
      +1.8000000e+002,      +6.3500000e-001,      +6.8890000e-001,
      +1.9000000e+002,      +4.6120000e-001,      +5.1190000e-001,
      +2.0000000e+002,      +3.4400000e-001,      +3.9040000e-001,
      +2.1000000e+002,      +2.6270000e-001,      +3.0640000e-001,
      +2.2000000e+002,      +2.0400000e-001,      +2.4510000e-001,
      +2.3000000e+002,      +1.5970000e-001,      +1.9910000e-001,
      +2.4000000e+002,      +1.2660000e-001,      +1.6380000e-001,
      +2.5000000e+002,      +1.0140000e-001,      +1.3600000e-001,
      +2.6000000e+002,      +8.1960000e-002,      +1.1400000e-001,
      +2.7000000e+002,      +6.6790000e-002,      +9.6270000e-002,
      +2.8000000e+002,      +5.4820000e-002,      +8.1950000e-002,
      +2.9000000e+002,      +4.5080000e-002,      +7.0250000e-002,
      +3.0000000e+002,      +3.7270000e-002,      +6.0530000e-002,
      +3.2000000e+002,      +2.5860000e-002,      +4.5540000e-002,
      +3.4000000e+002,      +1.8270000e-002,      +3.4770000e-002,
      +3.6000000e+002,      +1.3100000e-002,      +2.6890000e-002,
      +3.8000000e+002,      +9.5190000e-003,      +2.1010000e-002,
      +4.0000000e+002,      +6.9970000e-003,      +1.6580000e-002,
      +4.2000000e+002,      +5.1910000e-003,      +1.3170000e-002,
      +4.4000000e+002,      +3.8890000e-003,      +1.0560000e-002,
      +4.6000000e+002,      +2.9360000e-003,      +8.5180000e-003,
      +4.8000000e+002,      +2.2320000e-003,      +6.9150000e-003,
      +5.0000000e+002,      +1.7080000e-003,      +5.6450000e-003,
      +5.2000000e+002,      +1.3140000e-003,      +4.6310000e-003,
      +5.4000000e+002,      +1.0150000e-003,      +3.8180000e-003,
      +5.6000000e+002,      +7.8850000e-004,      +3.1610000e-003,
      +5.8000000e+002,      +6.1490000e-004,      +2.6270000e-003,
      +6.0000000e+002,      +4.8130000e-004,      +2.1920000e-003,
      +6.2000000e+002,      +3.7820000e-004,      +1.8350000e-003,
      +6.4000000e+002,      +2.9820000e-004,      +1.5410000e-003,
      +6.6000000e+002,      +2.3600000e-004,      +1.2980000e-003,
      +6.8000000e+002,      +1.8740000e-004,      +1.0960000e-003,
      +7.0000000e+002,      +1.4940000e-004,      +9.2810000e-004,
      +7.2000000e+002,      +1.1950000e-004,      +7.8780000e-004,
      +7.4000000e+002,      +9.5950000e-005,      +6.7020000e-004,
      +7.6000000e+002,      +7.7350000e-005,      +5.7140000e-004,
      +7.8000000e+002,      +6.2620000e-005,      +4.8810000e-004,
      +8.0000000e+002,      +5.0920000e-005,      +4.1790000e-004,
      +8.5000000e+002,      +3.1000000e-005,      +2.7000000e-004,
      +9.0000000e+002,      +1.8000000e-005,      +1.8000000e-004,
      +9.5000000e+002,      +1.1500000e-005,      +1.3000000e-004,
      +1.0000000e+003,      +7.6000000e-006,      +9.6000000e-005,
      +1.1000000e+003,      +3.8000000e-006,      +5.1000000e-005,
      +1.2000000e+003,      +2.1000000e-006,      +2.5000000e-005,
      +1.3000000e+003,      +1.3200000e-006,      +1.4000000e-005,
      +1.4000000e+003,      +9.4000000e-007,      +8.1000000e-006,
      +1.5000000e+003,      +7.0000000e-007,      +5.0000000e-006,
      +1.6000000e+003,      +5.2000000e-007,      +3.4500000e-006,
      +1.7000000e+003,      +4.1000000e-007,      +2.4000000e-006,
      +1.8000000e+003,      +3.3000000e-007,      +1.6200000e-006,
      +1.9000000e+003,      +2.6000000e-007,      +1.2000000e-006,
      +2.0000000e+003,      +2.1000000e-007,      +9.4000000e-007,
   },
   {   // F107 = +275.0000000
      +1.1000000e+002,      +7.8000000e+001,      +7.8000000e+001,
      +1.2000000e+002,      +2.4900000e+001,      +2.4900000e+001,
      +1.3000000e+002,      +7.7170000e+000,      +7.9840000e+000,
      +1.4000000e+002,      +3.7100000e+000,      +3.8410000e+000,
      +1.5000000e+002,      +2.1320000e+000,      +2.2150000e+000,
      +1.6000000e+002,      +1.3550000e+000,      +1.4200000e+000,
      +1.7000000e+002,      +9.2070000e-001,      +9.7750000e-001,
      +1.8000000e+002,      +6.5630000e-001,      +7.0520000e-001,
      +1.9000000e+002,      +4.8310000e-001,      +5.2800000e-001,
      +2.0000000e+002,      +3.6440000e-001,      +4.0640000e-001,
      +2.1000000e+002,      +2.8130000e-001,      +3.2010000e-001,
      +2.2000000e+002,      +2.2100000e-001,      +2.5740000e-001,
      +2.3000000e+002,      +1.7500000e-001,      +2.1020000e-001,
      +2.4000000e+002,      +1.4020000e-001,      +1.7370000e-001,
      +2.5000000e+002,      +1.1350000e-001,      +1.4510000e-001,
      +2.6000000e+002,      +9.2740000e-002,      +1.2220000e-001,
      +2.7000000e+002,      +7.6130000e-002,      +1.0380000e-001,
      +2.8000000e+002,      +6.2790000e-002,      +8.8910000e-002,
      +2.9000000e+002,      +5.2100000e-002,      +7.6560000e-002,
      +3.0000000e+002,      +4.3460000e-002,      +6.6250000e-002,
      +3.2000000e+002,      +3.0670000e-002,      +5.0260000e-002,
      +3.4000000e+002,      +2.2010000e-002,      +3.8690000e-002,
      +3.6000000e+002,      +1.6020000e-002,      +3.0150000e-002,
      +3.8000000e+002,      +1.1810000e-002,      +2.3750000e-002,
      +4.0000000e+002,      +8.8040000e-003,      +1.8870000e-002,
      +4.2000000e+002,      +6.6200000e-003,      +1.5110000e-002,
      +4.4000000e+002,      +5.0250000e-003,      +1.2190000e-002,
      +4.6000000e+002,      +3.8440000e-003,      +9.9060000e-003,
      +4.8000000e+002,      +2.9610000e-003,      +8.0910000e-003,
      +5.0000000e+002,      +2.2950000e-003,      +6.6540000e-003,
      +5.2000000e+002,      +1.7880000e-003,      +5.4970000e-003,
      +5.4000000e+002,      +1.4000000e-003,      +4.5610000e-003,
      +5.6000000e+002,      +1.1020000e-003,      +3.8010000e-003,
      +5.8000000e+002,      +8.7020000e-004,      +3.1800000e-003,
      +6.0000000e+002,      +6.9000000e-004,      +2.6710000e-003,
      +6.2000000e+002,      +5.4910000e-004,      +2.2510000e-003,
      +6.4000000e+002,      +4.3840000e-004,      +1.9020000e-003,
      +6.6000000e+002,      +3.5120000e-004,      +1.6130000e-003,
      +6.8000000e+002,      +2.8220000e-004,      +1.3710000e-003,
      +7.0000000e+002,      +2.2750000e-004,      +1.1690000e-003,
      +7.2000000e+002,      +1.8400000e-004,      +9.9840000e-004,
      +7.4000000e+002,      +1.4920000e-004,      +8.5500000e-004,
      +7.6000000e+002,      +1.2140000e-004,      +7.3380000e-004,
      +7.8000000e+002,      +9.9120000e-005,      +6.3100000e-004,
      +8.0000000e+002,      +8.1190000e-005,      +5.4370000e-004,
      +8.5000000e+002,      +4.4000000e-005,      +3.6500000e-004,
      +9.0000000e+002,      +2.7000000e-005,      +2.4000000e-004,
      +9.5000000e+002,      +1.7000000e-005,      +1.7500000e-004,
      +1.0000000e+003,      +1.1500000e-005,      +1.3000000e-004,
      +1.1000000e+003,      +5.3000000e-006,      +7.3000000e-005,
      +1.2000000e+003,      +2.6500000e-006,      +3.8000000e-005,
      +1.3000000e+003,      +1.5500000e-006,      +2.2000000e-005,
      +1.4000000e+003,      +1.0800000e-006,      +1.3500000e-005,
      +1.5000000e+003,      +8.0000000e-007,      +8.4000000e-006,
      +1.6000000e+003,      +5.8000000e-007,      +5.8000000e-006,
      +1.7000000e+003,      +4.6500000e-007,      +3.8000000e-006,
      +1.8000000e+003,      +3.6500000e-007,      +2.4000000e-006,
      +1.9000000e+003,      +2.9000000e-007,      +1.7000000e-006,
      +2.0000000e+003,      +2.3000000e-007,      +1.2500000e-006,
   } };

      // Constructor
   HarrisPriesterDrag::HarrisPriesterDrag()
      : workingF107(157.0)
   {
      init();
   }

   
   void HarrisPriesterDrag::init()
   {

      // load coefficient to easy accessing
      const int arrayF107[10] = {65, 75, 100, 125, 150, 175, 200, 225, 250, 275};
      for(int i = 0; i< 10; i++)
      {
         int f107 = arrayF107[i];

         Matrix<double> hpc(59,3,0.0);
         hpc= &hpcoef[i][0];

         hpCoeficient[f107] = hpc;
      }
      
   }  // End of method 'HarrisPriesterDrag::init()'


   Matrix<double> HarrisPriesterDrag::getDensityCoeficentsByF107(double f107)
   {
      map<int, Matrix<double> >::iterator itL,itH;
      
      itL = hpCoeficient.begin();

      map<int, Matrix<double> >::iterator it;
      for( it = hpCoeficient.begin(); it != hpCoeficient.end(); it++)
      {
         if( f107 >= it->first ) itL = it;
      }
      
      itH = itL;
      itH++;

      if( itH == hpCoeficient.end() )
      {
         itH = itL;
         itL--;
      }

      Matrix<double> out(59,3,0.0);
      for(int i=0; i<59; i++)
      {
         // H
         out[i][0] = itL->second[i][0];

         // m
         out[i][1] = itL->second[i][1] 
                   + (itH->second[i][1]-itL->second[i][1])
                      *(f107-itL->first)/(itH->first-itL->first);

         // M
         out[i][2] = itL->second[i][2] 
                   + (itH->second[i][2]-itL->second[i][2])
                      *(f107-itL->first)/(itH->first-itL->first);
      }

      return out;

   }  // End of method 'HarrisPriesterDrag::getDensityCoeficentsByF107()'


      // do some test here
   void HarrisPriesterDrag::test()
   {
      cout<<"testing HarrisPriesterDrag"<<endl;
      
      IERS::loadSTKFile("InputData\\EOP-v1.1.txt");
      ReferenceFrames::setJPLEphFile("InputData\\DE405\\jplde405");

      Vector<double> r(3),v(3);
      r(0)=-4453783.586;
      r(1)=-5038203.756;
      r(2)=-426384.456;

      v(0) =  3831.888;
      v(1) = -2887.221;
      v(2) = -6.018232;
      
      EarthBody body;
      UTCTime t;
      t.setMJD(53157.5);

      Spacecraft sc;
      sc.setDragArea(20.0);
      sc.setDragCoeff(2.2);
      sc.setDryMass(1000.0);

      Vector<double> rv(6,0.0);
      Vector<double> p(0,0);
      for(int i=0;i<3;i++)
      {
         rv(i) = r(i);
         rv(i+3) = v(i);
      }
      sc.initStateVector(rv);

      double den = computeDensity(t,body,r,v);
      doCompute(t,body,sc);
      
      Vector<double> accl = getAccel();
      
      double ax = accl(0);
      double ay = accl(1);
      double az = accl(2);

      int a = 0;
   }
   
   void HarrisPriesterDrag::updateF107(double f107)
   {
      workingF107 = f107;
      workingDens = getDensityCoeficentsByF107(workingF107);
   }

      /* Abstract class requires the subclass to compute the atmospheric density.
       * @param utc epoch in UTC
       * @param rb  EarthRef object.
       * @param r   Position vector.
       * @param v   Velocity vector
       * @return Atmospheric density in kg/m^3
       */
   double HarrisPriesterDrag::computeDensity( UTCTime utc, 
                                              EarthBody& rb, 
                                              Vector<double> r, 
                                              Vector<double> v)
   {
      double density = 0.0;
      
      // Get the J2000 to TOD transformation
      Matrix<double> N = ReferenceFrames::J2kToTODMatrix(utc);

      // Debuging
      /*
      double nn[3][3]={ {0.9999994803, -0.0009350126, -0.0004063480},
      {0.0009349995, 0.9999995624, -0.0000322758},
      {0.0004063780, 0.0000318958, 0.9999999169}};
      
      N = &nn[0][0];

      
      cout<<fixed<<setprecision(6);
      //cout<<workingDens<<endl;
      */

      // Transform r from J2000 to TOD
      Vector<double> r_tod = N * r;

      double rmag = norm(r_tod);
      
      //* Satellite true altitude
      Position pos(r_tod(0), r_tod(1), r_tod(2), Position::Cartesian);
      double alt = pos.getAltitude()/1000.0;      // km
      
      if ( alt >= upper_limit || alt <= lower_limit ) 
      { 
         //return 0.0;

         string msg = "HarrisPriesterDrag is good for 100.0 km t0 2000.0 km"
            + string("the altitude you try is ")
            + StringUtils::asString(alt) + " km!";
         
         Exception e(msg);
         
         //GPSTK_THROW(e);
      }

      Vector<double> r_Sun = ReferenceFrames::getJ2kPosition( utc.asTDB(),
                                                              SolarSystem::Sun);

      // get coefficients for this F107
      //updateF107(std::pow(149597870.0/norm(r_Sun),2)*dailyF107);
      updateF107(averageF107);

      // Debuging
      /*
      alt = 360.01323431364915;
      cout<<r_Sun<<endl;
      double rr[3] = {4.87202078904524E7,	1.318213179104784E8,	5.71498291483194E7};
      r_Sun = &rr[0];
      */

      double ra_Sun  = std::atan2( r_Sun(1), r_Sun(0));
      double dec_Sun = std::atan2( r_Sun(2), 
                                   std::sqrt( 
                                     std::pow(r_Sun(0),2) + std::pow(r_Sun(1),2)
                                     ) );

        //* Unit vector u towards the apex of the diurnal bulge
        //* in inertial geocentric coordinates
      double c_dec = std::cos(dec_Sun);
      

      Vector<double> u(3,0.0);      // Apex of diurnal bulge
      u(0) = c_dec * std::cos(ra_Sun + ra_lag);
      u(1) = c_dec * std::sin(ra_Sun + ra_lag);
      u(2) = std::sin(dec_Sun);

      //* Cosine of half angle between satellite position vector and
      //* apex of diurnal bulge
      double c_psi2 = (0.5 + 0.5 * dot(r,u)/norm(r));
      
      int ih = 0;

      // index search
      for(int i = 0; i < 59-1; i++)
      {
         if(alt >= workingDens[i][0] && alt < workingDens[i+1][0])
         {
            ih = i;
            break;
         }
      }
      
      // exponential density interpolation
      double h_min = (workingDens[ih][0] - workingDens[ih+1][0]) 
                   / std::log(workingDens[ih+1][1] / workingDens[ih][1]);

      double h_max = (workingDens[ih][0] - workingDens[ih+1][0]) 
                   / std::log(workingDens[ih+1][2] 
                   / workingDens[ih][2]);
      
      double d_min = workingDens[ih][1]*std::exp((workingDens[ih][0]-alt)/h_min);
      double d_max = workingDens[ih][2]*std::exp((workingDens[ih][0]-alt)/h_max);
      
      // Density computation
      //Use 2 for low inclination orbits and 6 for polar orbits
      // Reference Montenbruck P90
      double n_prm = 2.0;
      Vector<double> h = cross(r,v);
      double inc = std::acos(h(2) / norm(h));
      n_prm = 2.0 + inc * 8.0 / ASConstant::PI;


      density = d_min + (d_max - d_min) * std::pow(c_psi2,n_prm / 2.0);

      return density * 1.0e-9;   //[kg/m^3]

   }  // End of method 'HarrisPriesterDrag::computeDensity()'


}  // End of namespace 'gpstk'



